/*
 * BetaModel.cpp
 *
 * Copyright 2012 Pieter Agten
 *
 * This file is part of tltgen.
 *
 * tltgen is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * tltgen is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with tltgen.  If not, see <http://www.gnu.org/licenses/>.
 */

/**
 * @file   BetaModel.cpp
 * @author Pieter Agten (pieter.agten@gmail.com)
 * @date   4 sep 2012
 */

#include <cmath>
#include <stdexcept>
#include <iostream>
#include <list>
#include <boost/lexical_cast.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/make_shared.hpp>

#include "BetaModelThermistor.hpp"
#include "LeastSquaresSolver.hpp"
#include "ThermistorTable.hpp"
#include "ThermistorTableRow.hpp"

#include "BetaModel.hpp"

#define DEBUG false

using boost::shared_ptr;
using boost::make_shared;

namespace tltgen {

  BetaModel::BetaModel() { }

  BetaModel::~BetaModel() { }

  std::string BetaModel::getName() const
  {
    return "Beta model";
  }
  unsigned int BetaModel::getMinNbValues() const
  {
    return 2;
  }
  
  shared_ptr<Thermistor> BetaModel::generate(const shared_ptr<ThermistorTable>& thermistorTable, const LeastSquaresSolver& solver)
  {
    double params[3];
    solver.calculateParameters(this, thermistorTable.get(), params);

    Temperature T0 = params[0];
    Resistance R0 = params[1];
    double Beta = params[2];

    return make_shared<BetaModelThermistor>(T0, R0, Beta);
  }

  unsigned int BetaModel::getNbParameters() const
  {
    return 3;
  }

  std::string BetaModel::getParameterName(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument("Invalid parameter index specified.");
    }

    switch (pi) {
    case 0:
      return "T0";
    case 1:
      return "R0";
    case 2:
      return "Beta";
    default:
      throw std::runtime_error("Default switch case should not be reachable.");
    }
  }

  ThermistorModel::ParameterType BetaModel::getParameterType(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument("Invalid parameter index specified.");
    }

    switch (pi) {
    case 0:
      return ThermistorModel::TEMPERATURE;
    case 1:
      return ThermistorModel::RESISTANCE;
    case 2:
      return ThermistorModel::OTHER;
    default:
      throw std::runtime_error("Default switch case should not be reachable.");
    }
  }

  void BetaModel::getInitialGuess(double *guess) const
  {
    guess[0] = 298.15;
    guess[1] = 50000;
    guess[2] = 4000;
  }
  
  bool BetaModel::isParameterFixed(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument(boost::lexical_cast<std::string>(pi) + " is not a valid parameter index for this beta model generator.");
    }


    switch (pi) {
    case 0:
      // T0
      return true;
    case 1:
      // R0
      return false;
    case 2:
      // Beta
      return false;
    default:
      throw std::runtime_error("Default switch case should not be reachable.");
    }
  }

  bool BetaModel::hasParameterLowerBound(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument(boost::lexical_cast<std::string>(pi) + " is not a valid parameter index for this beta model generator.");
    }

    return true;
  }

  double BetaModel::getParameterLowerBound(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument(boost::lexical_cast<std::string>(pi) + " is not a valid parameter index for this beta model generator.");
    }

    switch (pi) {
    case 0:
      // Minimum T0
      return 0;
    case 1:
      // Minimum resistance value at T=T0
      return 1E-6;
    case 2:
      // Minimum beta value
      return 1;
    default:
      throw std::runtime_error("Default switch case should not be reachable.");
    }
  }

  bool BetaModel::hasParameterUpperBound(unsigned int pi) const
  {
    if (! isValidParameterIndex(pi)) {
      throw std::invalid_argument(boost::lexical_cast<std::string>(pi) + " is not a valid parameter index for this beta model generator.");
    }

    return false;
  }
  double BetaModel::getParameterUpperBound(unsigned int pi) const
  {
    throw std::runtime_error("This beta model generator has no upper bound on parameter " + boost::lexical_cast<std::string>(pi) + ".");
  }

  void BetaModel::applyModelFunction(const ThermistorTable* thermistorTable, const double *params, double *residues) const {
    Temperature t0 = params[0];
    Resistance r0 = params[1];
    double beta = params[2];

    BetaModelThermistor thermistor(t0, r0, beta);

#if DEBUG
    std::cout << std::endl << "[BMG] Trying parameters " << t0 << ", " << r0 << ", " << beta << "." << std::endl;
#endif
    unsigned int i = 0;
    for (std::list<ThermistorTableRow *>::const_iterator it = thermistorTable->begin(); it != thermistorTable->end(); ++it) {
      if ((*it)->isValid()) {
	double Ri = (*it)->getResistance();
	double Ti = thermistor.getTemperature(Ri);
	double residue = Ti - (*it)->getTemperature();
#if DEBUG
	std::cout << "[BMG] i=" << i << " R=" << Ri << " Ti=" << (*it)->getTemperature() << " Ti'=" << Ti << " Residue=" << residue << std::endl;
#endif
	residues[i] = residue;
	i += 1;
      }
    }
  }


  /**  
   * Jacobian matrix J(i,j) = dfi/dxj
   *   where fi = (Ti - ti)
   *         Ti = 1.0/(1.0/T0 + 1.0/B * ln(R/R0))
   *   the xj are the parameters R0 and B
   *
   *  dfi/dT0 = (Ti/T0)²     dfi/dR0 = Ti²/(B*R0)     dfi/dB = Ti²/B² * ln(R/R0)
   */
  void BetaModel::applyModelDerivative(const ThermistorTable* thermistorTable, const double *params, double *jacobian) const {
    Temperature t0 = params[0];
    Resistance r0 = params[1];
    double beta = params[2];

    BetaModelThermistor thermistor(t0, r0, beta);

    unsigned int i = 0;
    for (std::list<ThermistorTableRow *>::const_iterator it = thermistorTable->begin(); it != thermistorTable->end(); ++it) {
      if ((*it)->isValid()) {
	double Ri = (*it)->getResistance();
	double Ti = thermistor.getTemperature(Ri);
	double Ti2 = pow(Ti, 2);
	double B2 = pow(beta, 2);
      
	jacobian[i * 3 + 0] = pow(Ti/t0, 2);
	jacobian[i * 3 + 1] = Ti2/(beta * r0);
	jacobian[i * 3 + 2] = Ti2/B2 * log(Ri/r0);
	i += 1;
      }
    }
  }

}
